How does microbial community structure change with depth and oxygen concentration?

Alpha-diversity

Alpha-diversity and richness were calculated for the total community in R.

load("mothur_phyloseq.RData")

Samples were rarefied/normalized to 100,000 sequences per sample to facilitate comparisons between samples. A random seed was set to ensure reproducibility.

set.seed(4832)
m.norm = rarefy_even_depth(mothur, sample.size=100000)
## You set `rngseed` to FALSE. Make sure you've set & recorded
##  the random seed of your session for reproducibility.
## See `?set.seed`
## ...
## 614OTUs were removed because they are no longer 
## present in any sample after random subsampling
## ...

Rarefied counts were converted to relative abundance percentages.

m.perc = transform_sample_counts(m.norm, function(x) 100 * x/sum(x))
# Calculation
m.alpha = estimate_richness(m.norm, measures = c("Chao1", "Shannon"))
m.meta.alpha = full_join(rownames_to_column(m.alpha), rownames_to_column(data.frame(m.perc@sam_data)), by = "rowname")
m.meta.alpha

Example plots

m.meta.alpha %>% 
  
ggplot() +
  geom_point(aes(x=Depth_m, y=Shannon)) +
   geom_smooth(method='auto', aes(x=as.numeric(Depth_m), y=Shannon)) +
  labs(title="Alpha-diversity across depth", y="Shannon's diversity index", x="Depth (m)")
## `geom_smooth()` using method = 'loess'

m.meta.alpha %>% 

  ggplot() +
  geom_point(aes(x=O2_uM, y=Shannon)) +
   geom_smooth(method='auto', aes(x=as.numeric(O2_uM), y=Shannon)) +
  labs(title="Alpha-diversity across oxygen", y="Shannon's diversity index", x="O2 (uM)")
## `geom_smooth()` using method = 'loess'
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at -1.0833
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 33.437
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 1046.8
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used
## at -1.0833
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 33.437
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal
## condition number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other
## near singularities as well. 1046.8

m.meta.alpha %>% 

ggplot() +
  geom_point(aes(x=O2_uM, y=Shannon)) +
  labs(title="Alpha-diversity across oxygen", y="Shannon's diversity index", x="Oxygen (uM)")

m.meta.alpha %>% 
  mutate(O2_group = ifelse(O2_uM == 0, "anoxic", "oxic")) %>% 

# Use this in a plot
ggplot() +
  geom_boxplot(aes(x=O2_group, y=Shannon)) +
  labs(title="Example 3: Alpha-diversity by oxic/anoxic", y="Shannon's diversity index", x="Oxygen")

Taxa presence and abundance

Example plots using phyloseq for domain data. You should explore other taxonomy levels.

m.perc %>% 
plot_bar(fill="Domain") + 
  geom_bar(aes(fill=Domain), stat="identity") +
  labs(title="Example 4: Domains across samples") +
  labs(fill='Domain') +
  theme(plot.title = element_text(size=35, hjust=.5,vjust=.5,face="plain"),
        strip.text = element_text(size=12),
        axis.title.x = element_text(size=25),
        axis.title.y = element_text(size=25),
        axis.text.x = element_text(colour="grey20",size=12,angle=45,hjust=1,vjust=1),
        axis.text.y = element_text(colour="grey20",size=12),
        legend.text = element_text(size=12),
        legend.box.background = element_rect(),
        legend.box.margin = margin(6, 6, 6, 6),
        legend.position = "right")

m.perc %>% 
plot_bar() + 
  geom_bar(aes(fill=Domain), stat="identity") +
  facet_wrap(~Phylum, scales="free_y")+
  labs(title="Example 5: Phyla across samples")+
  labs(fill='Domain') +
  theme(plot.title = element_text(size=35, hjust=.5,vjust=.5,face="plain"),
        strip.text = element_text(size=5),
        axis.title.x = element_text(size=25),
        axis.title.y = element_text(size=25),
        axis.text.x = element_text(colour="grey20",size=9,angle=45,hjust=1,vjust=1),
        axis.text.y = element_text(colour="grey20",size=9),
        legend.text = element_text(size=12),
        legend.box.background = element_rect(),
        legend.box.margin = margin(6, 6, 6, 6),
        legend.position = "bottom")

Example plots using phyloseq for Family data.

# load subset
m.perc %>% 
plot_bar(fill="Family") + 
  geom_bar(aes(fill=Family), stat="identity") +
  labs(title="Example 6: Family across samples") +
  theme(plot.title = element_text(size=17.5, hjust=.5,vjust=.5,face="plain"),
        axis.title.x = element_text(size=12.5),
        axis.title.y = element_text(size=12.5),
        axis.text.x = element_text(colour="grey20",size=7.5,angle=45,hjust=1,vjust=1),
        axis.text.y = element_text(colour="grey20",size=7.5),
        legend.text = element_text(size=6),
        legend.box.background = element_rect(),
        legend.box.margin = margin(6, 6, 6, 6),
        legend.position = "bottom")

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
plot_bar(fill="Genus") + 
  geom_bar(aes(fill=Genus), stat="identity") +
  labs(title="Example 7: Genus of family Oceanospirillaceae across samples") +
  labs(fill="Genus") +
  theme(plot.title = element_text(size=17.5, hjust=.5,vjust=.5,face="plain"),
        axis.title.x = element_text(size=12.5),
        axis.title.y = element_text(size=12.5),
        axis.text.x = element_text(colour="grey20",size=7.5,angle=45,hjust=1,vjust=1),
        axis.text.y = element_text(colour="grey20",size=7.5),
        legend.text = element_text(size=6),
        legend.box.background = element_rect(),
        legend.box.margin = margin(6, 6, 6, 6),
        legend.position = "bottom")

Examples outside of phyloseq

m.perc %>%
  tax_glom(taxrank = 'Domain') %>%
  psmelt() %>% 

ggplot() +
  geom_boxplot(aes(x=Domain, y=Abundance)) +
  coord_flip() +
  labs(title="Example 7: Domain boxplots")+
  theme(plot.title = element_text(size=17.5, hjust=.5,vjust=.5,face="plain"),
      axis.title.x = element_text(size=12.5),
      axis.title.y = element_text(size=12.5)
  )

m.perc %>% 
  tax_glom(taxrank = 'Family') %>%
  psmelt() %>% 

ggplot() +
  geom_boxplot(aes(x=Family, y=Abundance)) +
  coord_flip() +
  labs(title="Example 7: Family boxplots")+
  theme(plot.title = element_text(size=17.5, hjust=.5,vjust=.5,face="plain"),
    axis.title.x = element_text(size=12.5),
    axis.title.y = element_text(size=12.5)
  )

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
  psmelt() %>% 

ggplot() +
  geom_boxplot(aes(x=Genus, y=Abundance)) +
  coord_flip() +
  labs(title="Example 7: Domain boxplots")+
  theme(plot.title = element_text(size=17.5, hjust=.5,vjust=.5,face="plain"),
    axis.title.x = element_text(size=12.5),
    axis.title.y = element_text(size=12.5)
  )

Does your taxon of interest significantly differ in abundance with depth and/or oxygen concentration?

Using the magrittr package, we can pipe our tidyverse modified data into linear models and other statiscal tests.

Depth Linear model

m.norm %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
  tax_glom(taxrank = 'Family') %>%
  psmelt() %>%

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       6       5       7 
##  -71.71  379.84   44.85 -193.39  -19.87 -293.63  153.90 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1892.219    227.762   8.308 0.000413 ***
## Depth_m      -10.051      1.656  -6.070 0.001753 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 244.6 on 5 degrees of freedom
## Multiple R-squared:  0.8805, Adjusted R-squared:  0.8566 
## F-statistic: 36.84 on 1 and 5 DF,  p-value: 0.001753

Plot to go along with linear model above.

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
  psmelt() %>% 
  group_by(Sample) %>% 
  summarize(Abundance_sum=sum(Abundance), Depth_m=mean(Depth_m)) %>% 
  
ggplot() +
  geom_point(aes(x=Depth_m, y=Abundance_sum)) +
  geom_smooth(method='lm', aes(x=as.numeric(Depth_m), y=Abundance_sum)) +
  labs(title="Example 9: Abundance of Family Oceanospirillaceae across depth")

Oxygen Linear model

m.norm %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
  tax_glom(taxrank = 'Family') %>%
  psmelt() %>%

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##      1      2      3      4      6      5      7 
## -135.6  680.3  184.5 -120.0 -102.7 -225.7 -280.7 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  316.742    162.027   1.955    0.108  
## O2_uM          7.102      1.920   3.699    0.014 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 366 on 5 degrees of freedom
## Multiple R-squared:  0.7324, Adjusted R-squared:  0.6789 
## F-statistic: 13.68 on 1 and 5 DF,  p-value: 0.01401

Plot to go along with linear model above.

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
  psmelt() %>% 
  group_by(Sample) %>% 
  summarize(Abundance_sum=sum(Abundance), O2_uM=mean(O2_uM)) %>% 
  
ggplot() +
  geom_point(aes(x=O2_uM, y=Abundance_sum)) +
  geom_smooth(method='lm', aes(x=as.numeric(O2_uM), y=Abundance_sum)) +
  labs(title="Example 10: Abundance of Family Oceanospirillaceaes across oxygen")

Within your taxon, what is the richness (number of OTUs/ASVs)?

Across all samples, there are 25 OTUs (i.e. taxa) within family Oceanospirillaceae

m.norm %>% 
  subset_taxa(Family=="Oceanospirillaceae") 
## phyloseq-class experiment-level object
## otu_table()   OTU Table:         [ 25 taxa and 7 samples ]
## sample_data() Sample Data:       [ 7 samples by 22 sample variables ]
## tax_table()   Taxonomy Table:    [ 25 taxa by 7 taxonomic ranks ]

Counting OTUs within family Oceanospirillaceae within each sample…

m.norm %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>%
  estimate_richness(measures = c("Observed"))

Do the abundances of OTUs/ASVs within your taxon of interest change significantly with depth and/or oxygen concentration?

Depth General linear model for each OTU

#OTU0065
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0065") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary() 
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  298.09 -331.31 -215.40 -110.47  -32.54   72.40  319.24 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept) 1119.865    251.744   4.448  0.00671 **
## Depth_m       -7.196      1.830  -3.932  0.01105 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 270.3 on 5 degrees of freedom
## Multiple R-squared:  0.7556, Adjusted R-squared:  0.7067 
## F-statistic: 15.46 on 1 and 5 DF,  p-value: 0.01105
#OTU0084
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0084") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary() 
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
## 203.974 -97.115  49.328   3.343 -77.643 -63.628 -18.260 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) 311.7915   105.2621   2.962   0.0314 *
## Depth_m      -1.4677     0.7652  -1.918   0.1132  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 113 on 5 degrees of freedom
## Multiple R-squared:  0.4239, Adjusted R-squared:  0.3086 
## F-statistic: 3.678 on 1 and 5 DF,  p-value: 0.1132
#OTU0090
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0090") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  159.66   98.37   51.97  -25.10  -75.56 -106.55 -102.78 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept) 114.90573  107.08535   1.073    0.332
## Depth_m      -0.03561    0.77849  -0.046    0.965
## 
## Residual standard error: 115 on 5 degrees of freedom
## Multiple R-squared:  0.0004184,  Adjusted R-squared:  -0.1995 
## F-statistic: 0.002093 on 1 and 5 DF,  p-value: 0.9653
#OTU0104
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0104") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  200.76   92.39  -32.39  -37.95  -49.17  -60.10 -113.56 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) 121.3728   110.3807   1.100    0.322
## Depth_m      -0.2814     0.8024  -0.351    0.740
## 
## Residual standard error: 118.5 on 5 degrees of freedom
## Multiple R-squared:  0.024,  Adjusted R-squared:  -0.1712 
## F-statistic: 0.123 on 1 and 5 DF,  p-value: 0.7401
#OTU0117
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0117") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##      1      2      3      4      5      6      7 
## 171.73  38.05  -7.46 -84.20 -44.97 -42.48 -30.68 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) 104.8550    86.2072   1.216    0.278
## Depth_m      -0.3659     0.6267  -0.584    0.585
## 
## Residual standard error: 92.58 on 5 degrees of freedom
## Multiple R-squared:  0.06382,    Adjusted R-squared:  -0.1234 
## F-statistic: 0.3409 on 1 and 5 DF,  p-value: 0.5847
#OTU0327
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0327") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  15.064 -16.764  15.762  -6.030  -1.001 -11.059   4.028 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept) 51.28903   12.71492   4.034  0.00998 **
## Depth_m     -0.33525    0.09244  -3.627  0.01511 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 13.65 on 5 degrees of freedom
## Multiple R-squared:  0.7246, Adjusted R-squared:  0.6695 
## F-statistic: 13.15 on 1 and 5 DF,  p-value: 0.01511
#OTU0418
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0418") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##      1      2      3      4      5      6      7 
##  9.921  8.646 -1.460 -1.872 -3.979 -4.666 -6.591 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept)  6.72831    6.71325   1.002    0.362
## Depth_m     -0.01375    0.04880  -0.282    0.789
## 
## Residual standard error: 7.209 on 5 degrees of freedom
## Multiple R-squared:  0.01562,    Adjusted R-squared:  -0.1813 
## F-statistic: 0.07935 on 1 and 5 DF,  p-value: 0.7895
#OTU0511
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0511") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  6.9555 -6.4543 -5.9899 -0.7931  7.8681 -3.3915  1.8052 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept) 26.77676    5.93643   4.511  0.00634 **
## Depth_m     -0.17322    0.04316  -4.014  0.01018 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.375 on 5 degrees of freedom
## Multiple R-squared:  0.7632, Adjusted R-squared:  0.7158 
## F-statistic: 16.11 on 1 and 5 DF,  p-value: 0.01018
#OTU0675
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0675") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  4.1476  4.6213 -1.6723 -0.4507 -4.5227 -2.8939  0.7709 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) 12.66678    3.50856   3.610   0.0154 *
## Depth_m     -0.08144    0.02551  -3.193   0.0242 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.768 on 5 degrees of freedom
## Multiple R-squared:  0.6709, Adjusted R-squared:  0.6051 
## F-statistic: 10.19 on 1 and 5 DF,  p-value: 0.02419
#OTU0857
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0857") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  9.0638 -4.1489 -2.6170 -3.3830  3.7660 -1.8511 -0.8298 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.27660    4.85570  -0.881    0.419
## Depth_m      0.05106    0.03530   1.447    0.208
## 
## Residual standard error: 5.215 on 5 degrees of freedom
## Multiple R-squared:  0.295,  Adjusted R-squared:  0.154 
## F-statistic: 2.093 on 1 and 5 DF,  p-value: 0.2077
#OTU0952
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0952") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  1.21768 -0.07038  0.32733  1.25532 -0.46809 -1.39607 -0.86579 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  4.04746    1.02755   3.939   0.0110 *
## Depth_m     -0.02651    0.00747  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.103 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU0979
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0979") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  1.9483 -0.1126  2.0085  0.5237 -0.7489 -2.2337 -1.3853 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  6.47594    1.64407   3.939   0.0110 *
## Depth_m     -0.04242    0.01195  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.766 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU1077
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1077") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  1.46121  0.39280 -0.08445 -0.56170 -1.03895  1.50638 -1.67529 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  4.856956   1.233055   3.939   0.0110 *
## Depth_m     -0.031817   0.008964  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.324 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU1349
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1349") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.73061 -0.28085 -0.04223  0.75319 -0.51948  0.19640 -0.83764 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  2.428478   0.616528   3.939   0.0110 *
## Depth_m     -0.015908   0.004482  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6621 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU1516
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1516") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.9741 -0.6926  1.0043 -0.3745 -0.0563 -1.1169  0.2619 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  3.237971   0.822037   3.939   0.0110 *
## Depth_m     -0.021211   0.005976  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.8828 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU1678
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1678") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.9617  1.4645 -0.6298 -0.8183  0.1558 -0.7241 -0.4098 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.218658   0.917334  -0.238    0.821
## Depth_m      0.006285   0.006669   0.942    0.389
## 
## Residual standard error: 0.9851 on 5 degrees of freedom
## Multiple R-squared:  0.1508, Adjusted R-squared:  -0.019 
## F-statistic: 0.8881 on 1 and 5 DF,  p-value: 0.3893
#OTU1685
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1685") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.73061 -0.04223  0.75319 -0.51948 -0.28085 -0.83764  0.19640 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  2.428478   0.616528   3.939   0.0110 *
## Depth_m     -0.015908   0.004482  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6621 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU1730
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1730") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  2.2511  2.2592  0.1666 -1.3961 -1.0936 -1.2448 -0.9424 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.26743    1.65437  -0.162    0.878
## Depth_m      0.01008    0.01203   0.838    0.440
## 
## Residual standard error: 1.777 on 5 degrees of freedom
## Multiple R-squared:  0.1232, Adjusted R-squared:  -0.05213 
## F-statistic: 0.7027 on 1 and 5 DF,  p-value: 0.4401
#OTU1992
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1992") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.24354 -0.17316 -0.01408  0.25106  0.06547 -0.09362 -0.27921 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  0.809493   0.205509   3.939   0.0110 *
## Depth_m     -0.005303   0.001494  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2207 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU2430
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu2430") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.82684 -0.11424 -0.14959 -0.27921 -0.09656 -0.05532 -0.13191 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.290998   0.378670   0.768    0.477
## Depth_m     -0.001178   0.002753  -0.428    0.686
## 
## Residual standard error: 0.4067 on 5 degrees of freedom
## Multiple R-squared:  0.03535,    Adjusted R-squared:  -0.1576 
## F-statistic: 0.1832 on 1 and 5 DF,  p-value: 0.6864
#OTU3469
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3469") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.82684 -0.13191 -0.09656 -0.14959 -0.27921 -0.11424 -0.05532 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.290998   0.378670   0.768    0.477
## Depth_m     -0.001178   0.002753  -0.428    0.686
## 
## Residual standard error: 0.4067 on 5 degrees of freedom
## Multiple R-squared:  0.03535,    Adjusted R-squared:  -0.1576 
## F-statistic: 0.1832 on 1 and 5 DF,  p-value: 0.6864
#OTU3676
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3676") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  1.2085 -0.4511 -0.5532 -0.3489 -0.1106 -0.2468  0.5021 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.570213   0.647427  -0.881    0.419
## Depth_m      0.006809   0.004707   1.447    0.208
## 
## Residual standard error: 0.6953 on 5 degrees of freedom
## Multiple R-squared:  0.295,  Adjusted R-squared:  0.154 
## F-statistic: 2.093 on 1 and 5 DF,  p-value: 0.2077
#OTU3677
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3677") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.73061 -0.28085  0.19640 -0.51948  0.75319 -0.04223 -0.83764 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  2.428478   0.616528   3.939   0.0110 *
## Depth_m     -0.015908   0.004482  -3.549   0.0164 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6621 on 5 degrees of freedom
## Multiple R-squared:  0.7159, Adjusted R-squared:  0.6591 
## F-statistic:  12.6 on 1 and 5 DF,  p-value: 0.0164
#OTU3678
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3678") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.82684 -0.13191 -0.05532 -0.11424 -0.14959 -0.09656 -0.27921 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.290998   0.378670   0.768    0.477
## Depth_m     -0.001178   0.002753  -0.428    0.686
## 
## Residual standard error: 0.4067 on 5 degrees of freedom
## Multiple R-squared:  0.03535,    Adjusted R-squared:  -0.1576 
## F-statistic: 0.1832 on 1 and 5 DF,  p-value: 0.6864
#OTU3783
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3783") %>% 

  lm(Abundance ~ Depth_m, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ Depth_m, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.60426 -0.05532 -0.27660 -0.12340 -0.17447  0.25106 -0.22553 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.285106   0.323713  -0.881    0.419
## Depth_m      0.003404   0.002353   1.447    0.208
## 
## Residual standard error: 0.3476 on 5 degrees of freedom
## Multiple R-squared:  0.295,  Adjusted R-squared:  0.154 
## F-statistic: 2.093 on 1 and 5 DF,  p-value: 0.2077

Repeated for all OTUs within domain and then corrected for multiple comparisons.

p.adjust(c(0.01105, 0.1132, 0.9653, 0.7401, 0.5847, 0.01511, 0.7895, 0.01018, 0.02419, 0.2077, 0.0164, 0.0164, 0.0164, 0.0164, 0.0164, 0.3893, 0.0164, 0.4401, 0.0164, 0.6864, 0.6864, 0.2077, 0.0164, 0.6864, 0.2077), method="fdr")
##  [1] 0.03727273 0.21769231 0.96530000 0.80445652 0.76934211 0.03727273
##  [7] 0.82239583 0.03727273 0.05039583 0.32453125 0.03727273 0.03727273
## [13] 0.03727273 0.03727273 0.03727273 0.57250000 0.03727273 0.61125000
## [19] 0.03727273 0.78000000 0.78000000 0.32453125 0.03727273 0.78000000
## [25] 0.32453125

Plots to go along with these tests.

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
  psmelt() %>% 
  
ggplot() +
  geom_point(aes(x=Depth_m, y=Abundance)) +
  geom_smooth(method='lm', aes(x=Depth_m, y=Abundance)) +
  facet_wrap(~OTU, scales="free_y") +
  labs(title="Example 11: Abundance of OTUs within Family Oceanospirillaceae across depth")

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>%
  psmelt() %>% 
  
ggplot() +
  geom_point(aes(x=Sample, y=OTU, size=Abundance, color=OTU)) + 
  scale_size_continuous(range = c(0,5)) +
  labs(title="Abundance of OTUs within Family Oceanospirillaceae across depth")+
  theme(plot.title = element_text(size=15, hjust=.5,vjust=.5,face="plain"),
        axis.text.x = element_text(colour="grey20",size=10,angle=45,hjust=.5,vjust=.5,face="plain"),
        legend.box.background = element_rect(),
        legend.box.margin = margin(6, 6, 6, 6),
        legend.position = "right")

Oxygen

#OTU0065
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0065") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##   36.88 -108.84 -101.01  -28.61   70.86   67.86   62.86 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -62.8605    38.4906  -1.633    0.163    
## O2_uM         6.3322     0.4561  13.884 3.48e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 86.95 on 5 degrees of freedom
## Multiple R-squared:  0.9747, Adjusted R-squared:  0.9697 
## F-statistic: 192.8 on 1 and 5 DF,  p-value: 3.483e-05
#OTU0084
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0084") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
## 246.042 -53.699  66.182   6.897 -81.141 -89.141 -95.141 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept)  95.1407    59.8576   1.589    0.173
## O2_uM         0.7318     0.7093   1.032    0.349
## 
## Residual standard error: 135.2 on 5 degrees of freedom
## Multiple R-squared:  0.1755, Adjusted R-squared:  0.01064 
## F-statistic: 1.065 on 1 and 5 DF,  p-value: 0.3495
#OTU0090
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0090") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  158.39   94.30   34.35  -34.10  -92.65  -38.65 -121.65 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) 126.6462    48.9713   2.586   0.0491 *
## O2_uM        -0.3692     0.5803  -0.636   0.5526  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 110.6 on 5 degrees of freedom
## Multiple R-squared:  0.0749, Adjusted R-squared:  -0.1101 
## F-statistic: 0.4048 on 1 and 5 DF,  p-value: 0.5526
#OTU0104
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0104") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##      1      2      3      4      5      6      7 
## 206.81  91.67 -39.73 -57.84 -64.84 -89.84 -46.22 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept)  94.8434    52.5659   1.804    0.131
## O2_uM        -0.2013     0.6229  -0.323    0.760
## 
## Residual standard error: 118.8 on 5 degrees of freedom
## Multiple R-squared:  0.02047,    Adjusted R-squared:  -0.1754 
## F-statistic: 0.1045 on 1 and 5 DF,  p-value: 0.7596
#OTU0117
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0117") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##      1      2      3      4      5      6      7 
## 180.96  39.79 -11.57 -36.60 -55.19 -58.19 -59.19 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) 60.19451   42.33796   1.422    0.214
## O2_uM       -0.03045    0.50168  -0.061    0.954
## 
## Residual standard error: 95.65 on 5 degrees of freedom
## Multiple R-squared:  0.0007361,  Adjusted R-squared:  -0.1991 
## F-statistic: 0.003683 on 1 and 5 DF,  p-value: 0.954
#OTU0327
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0327") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##      1      2      3      4      5      6      7 
##  2.166 -6.373  4.001 -2.117  4.001 -5.680  4.001 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -4.00141    2.25787  -1.772    0.137    
## O2_uM        0.29924    0.02675  11.185 9.97e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5.101 on 5 degrees of freedom
## Multiple R-squared:  0.9616, Adjusted R-squared:  0.9539 
## F-statistic: 125.1 on 1 and 5 DF,  p-value: 9.97e-05
#OTU0418
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0418") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##      1      2      3      4      5      6      7 
##  9.854  8.925 -2.554 -2.296 -5.554 -5.554 -2.821 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept)  5.55418    3.18102   1.746    0.141
## O2_uM       -0.01262    0.03769  -0.335    0.751
## 
## Residual standard error: 7.186 on 5 degrees of freedom
## Multiple R-squared:  0.02192,    Adjusted R-squared:  -0.1737 
## F-statistic: 0.112 on 1 and 5 DF,  p-value: 0.7514
#OTU0511
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0511") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.8143 -1.1036 -3.2459  1.6586  1.6586 -1.4407  1.6586 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.65864    0.94412  -1.757    0.139    
## O2_uM        0.15159    0.01119  13.550 3.92e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.133 on 5 degrees of freedom
## Multiple R-squared:  0.9735, Adjusted R-squared:  0.9682 
## F-statistic: 183.6 on 1 and 5 DF,  p-value: 3.922e-05
#OTU0675
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0675") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.6452  1.8584 -0.6716  0.8584 -1.9861 -1.5627  0.8584 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.858387   0.690901  -1.242 0.269173    
## O2_uM        0.074830   0.008187   9.140 0.000263 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.561 on 5 degrees of freedom
## Multiple R-squared:  0.9435, Adjusted R-squared:  0.9322 
## F-statistic: 83.55 on 1 and 5 DF,  p-value: 0.0002626
#OTU0857
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0857") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
## 12.0609 -2.9391 -2.5685 -2.9391  0.9884 -2.3526 -2.2501 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept)  2.93910    2.66264   1.104     0.32
## O2_uM       -0.01813    0.03155  -0.575     0.59
## 
## Residual standard error: 6.015 on 5 degrees of freedom
## Multiple R-squared:  0.06193,    Adjusted R-squared:  -0.1257 
## F-statistic: 0.3301 on 1 and 5 DF,  p-value: 0.5905
#OTU0952
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0952") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.1810  0.3295  0.3295  0.3295 -0.1564 -0.5738 -0.4393 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.329477   0.188332  -1.749 0.140622    
## O2_uM        0.023762   0.002232  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4255 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU0979
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu0979") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.2896  0.5272  0.5272  0.5272 -0.2502 -0.9180 -0.7029 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.527163   0.301331  -1.749 0.140622    
## O2_uM        0.038019   0.003571  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6807 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU1077
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1077") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.2172  0.3954  0.3954 -0.1876 -0.5272  0.3954 -0.6885 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.395373   0.225998  -1.749 0.140622    
## O2_uM        0.028514   0.002678  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5106 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU1349
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1349") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.10861 -0.09382  0.19769  0.19769 -0.26359  0.19769 -0.34426 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.197686   0.112999  -1.749 0.140622    
## O2_uM        0.014257   0.001339  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2553 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU1516
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1516") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.1448 -0.3515  0.2636 -0.1251  0.2636 -0.4590  0.2636 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.263582   0.150665  -1.749 0.140622    
## O2_uM        0.019010   0.001785  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3404 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU1678
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1678") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  1.29444  1.39324 -0.64313 -0.70556 -0.04394 -0.70556 -0.58949 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.705562   0.458965   1.537    0.185
## O2_uM       -0.003054   0.005438  -0.561    0.599
## 
## Residual standard error: 1.037 on 5 degrees of freedom
## Multiple R-squared:  0.05931,    Adjusted R-squared:  -0.1288 
## F-statistic: 0.3153 on 1 and 5 DF,  p-value: 0.5987
#OTU1685
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1685") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.10861  0.19769  0.19769 -0.26359 -0.09382 -0.34426  0.19769 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.197686   0.112999  -1.749 0.140622    
## O2_uM        0.014257   0.001339  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2553 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU1730
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1730") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  2.76623  1.96853 -0.08068 -1.23377 -1.12496 -1.23377 -1.06158 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept)  1.233770   0.815500   1.513    0.191
## O2_uM       -0.005322   0.009663  -0.551    0.606
## 
## Residual standard error: 1.842 on 5 degrees of freedom
## Multiple R-squared:  0.05719,    Adjusted R-squared:  -0.1314 
## F-statistic: 0.3033 on 1 and 5 DF,  p-value: 0.6055
#OTU1992
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu1992") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.03620 -0.08786  0.06590  0.06590  0.06590 -0.03127 -0.11475 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0658954  0.0376663  -1.749 0.140622    
## O2_uM        0.0047524  0.0004463  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.08509 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU2430
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu2430") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.8562 -0.1500 -0.1447 -0.1148 -0.1500 -0.1500 -0.1467 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.1500035  0.1831720   0.819    0.450
## O2_uM       -0.0001627  0.0021705  -0.075    0.943
## 
## Residual standard error: 0.4138 on 5 degrees of freedom
## Multiple R-squared:  0.001122,   Adjusted R-squared:  -0.1987 
## F-statistic: 0.005619 on 1 and 5 DF,  p-value: 0.9432
#OTU3469
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3469") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.8562 -0.1467 -0.1500 -0.1447 -0.1148 -0.1500 -0.1500 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.1500035  0.1831720   0.819    0.450
## O2_uM       -0.0001627  0.0021705  -0.075    0.943
## 
## Residual standard error: 0.4138 on 5 degrees of freedom
## Multiple R-squared:  0.001122,   Adjusted R-squared:  -0.1987 
## F-statistic: 0.005619 on 1 and 5 DF,  p-value: 0.9432
#OTU3676
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3676") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  1.6081 -0.3919 -0.3919 -0.3425 -0.3000 -0.3137  0.1318 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.391880   0.355018   1.104     0.32
## O2_uM       -0.002417   0.004207  -0.575     0.59
## 
## Residual standard error: 0.802 on 5 degrees of freedom
## Multiple R-squared:  0.06193,    Adjusted R-squared:  -0.1257 
## F-statistic: 0.3301 on 1 and 5 DF,  p-value: 0.5905
#OTU3677
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3677") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##        1        2        3        4        5        6        7 
##  0.10861 -0.09382  0.19769 -0.26359  0.19769  0.19769 -0.34426 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.197686   0.112999  -1.749 0.140622    
## O2_uM        0.014257   0.001339  10.648 0.000126 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2553 on 5 degrees of freedom
## Multiple R-squared:  0.9578, Adjusted R-squared:  0.9493 
## F-statistic: 113.4 on 1 and 5 DF,  p-value: 0.0001264
#OTU3678
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3678") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.8562 -0.1467 -0.1500 -0.1500 -0.1447 -0.1500 -0.1148 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.1500035  0.1831720   0.819    0.450
## O2_uM       -0.0001627  0.0021705  -0.075    0.943
## 
## Residual standard error: 0.4138 on 5 degrees of freedom
## Multiple R-squared:  0.001122,   Adjusted R-squared:  -0.1987 
## F-statistic: 0.005619 on 1 and 5 DF,  p-value: 0.9432
#OTU3783
m.norm %>% 
  psmelt() %>% 
  filter(OTU=="Otu3783") %>% 

  lm(Abundance ~ O2_uM, .) %>% 
  summary()
## 
## Call:
## lm(formula = Abundance ~ O2_uM, data = .)
## 
## Residuals:
##       1       2       3       4       5       6       7 
##  0.8041 -0.1500 -0.1959 -0.1568 -0.1712  0.0659 -0.1959 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.195940   0.177509   1.104     0.32
## O2_uM       -0.001208   0.002103  -0.575     0.59
## 
## Residual standard error: 0.401 on 5 degrees of freedom
## Multiple R-squared:  0.06193,    Adjusted R-squared:  -0.1257 
## F-statistic: 0.3301 on 1 and 5 DF,  p-value: 0.5905
p.adjust(c(0.00003483, 0.3495, 0.5526, 0.7596, 0.954, 0.0000997, 0.7514, 0.00003922, 0.0002626, 0.5905, 0.0001264, 0.0001264, 0.0001264, 0.0001264, 0.0001264, 0.5987, 0.0001264, 0.6055, 0.0001264, 0.9432, 0.9432, 0.5905, 0.0164, 0.9432, 0.5905), method="fdr")
##  [1] 0.0003160000 0.6721153846 0.7967105263 0.9042857143 0.9540000000
##  [6] 0.0003160000 0.9042857143 0.0003160000 0.0005968182 0.7967105263
## [11] 0.0003160000 0.0003160000 0.0003160000 0.0003160000 0.0003160000
## [16] 0.7967105263 0.0003160000 0.7967105263 0.0003160000 0.9540000000
## [21] 0.9540000000 0.7967105263 0.0341666667 0.9540000000 0.7967105263

Plots to go along with these tests.

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>% 
  psmelt() %>% 
  
ggplot() +
  geom_point(aes(x=O2_uM, y=Abundance)) +
  geom_smooth(method='lm', aes(x=O2_uM, y=Abundance)) +
  facet_wrap(~OTU, scales="free_y") +
  labs(title="Example 11: Abundance of OTUs within Family Oceanospirillaceae across oxygen")

m.perc %>% 
  subset_taxa(Family=="Oceanospirillaceae") %>%
  psmelt() %>% 
  
ggplot() +
  geom_point(aes(x=O2_uM, y=OTU, size=Abundance, color=OTU)) + 
  scale_size_continuous(range = c(0,5)) +
  labs(title="Example 12: Abundance of OTUs within Family Oceanospirillaceae across oxygen")+
    theme(plot.title = element_text(size=15, hjust=.5,vjust=.5,face="plain"),
        axis.text.x = element_text(colour="grey20",size=10,angle=0,hjust=.5,vjust=.5,face="plain"),
        legend.box.background = element_rect(),
        legend.box.margin = margin(6, 6, 6, 6),
        legend.position = "right")